Computational experience with DICOPT solving MINLP problems in process systems engineering

This paper discusses an implementation of the Outer-Approximation/Equality-Relaxation algorithm for solving MINLP problems that arise in process systems engineering. The computer code DICOPT has been developed utilizing state-of-the-art optimization tools and a powerful modelling language. Computational experience in solving fifteen MINLP problems with DICOPT is reported. Applications include design of batch processes, structural flowsheet optimization, planning, flexibility, and reliability problems. Results show that DICOPT provides a very efficient tool for solving MINLP problems. INTRODUCTION Many optimization models for design and operation of process systems require the treatment of both continuous and discrete (mainly 0-1) variables. The presence of nonlinearities in the models gives rise to mixed-integer nonlinear programming (MINLP) problems, which in the past have proved to be very expensive and difficult to solve. It is shown in this paper, however, that this situation has undergone a considerable change and that reliable and efficient algorithms now exist The development of the Outer-Approximation (OA) algorithm by Duran and Grossmann (1986a), its extension with the Equality-Relaxation (OA/ER) strategy by Kocis and Grossmann (1987a), coupled with advances in nonlinear program (NLP) solvers (e.g. MINOS, Murtagh and Saunders, 1985) and mixed-integer linear program (MILP) solvers (e.g. MPSX, IBM, 1979), has led to the capability of solving MINLP problems quite efficiently. These algorithms and techniques have been implemented in the computer program DICOPT (Discrete Continuous OPTimizer) that features as interface the modelling language GAMS (General Algebraic Modelling System, Kendrick and Meeraus, 1985) which greatly facilitates the algebraic formulation of MINLP problems. Computational results are reported on fifteen MINLP problems that have been solved with DICOPT. These include applications to design of batch processes, structural flowsheet optimization, planning, flexibility, and reliability design. Computer times on an IBM-3090/600 are typically of the order of only 2-100 seconds for problems with up to 60 0-1 variables and 410 variables and constraints. This paper also briefly discusses a number of points that arise from this study and that seem to be rather unique to MINLP problems. The points include the effect that different modelling schemes have not only on the algorithm's efficiency but also on the global optimality of the solution, the role of convexity in the quality of the predicted bounds, and aspects of computer implementation that have to do with the use of different optimizers. Future directions are also briefly discussed. OUTER-APPROXIMATION/EQUALITY-RELAXATION ALGORITHM The OA/ER algorithm addresses MINLP problems of the general form: z =min cy+f(x) s.t. h(x)=0 g(x)<0 Ax=a (MINLP) By+Cx(x,y), can be treated through reformulation to yield a problem which is linear in y. This can be done simply by defining new continuous variables x = y, so that \|/(y) and <(>(x,y) can be formulated as nonlinear functions of continuous variables V(x) and (Ĵ XjX). The additional equations x = y are then linear in y. This, however, is equivalent to treating y(y) and <|>(x,y) directly as nonlinear functions which are linearized in x,y at the NLP subproblem solution points x,y, k=l,2,...K, to define the master problem. It should be noted that the OA/ER is not the only method available for solving MINLP problems. Branch and bound techniques can be applied to MINLP problems as well as the Generalized Benders Decomposition (GBD) (Benders, 1962, Geoffrion, 1972). The latter is an iterative bounding algorithm which is similar to OA/ER except that its master problem is based on a dual representation of the NLP subproblems, and generally predicts weak lower bounds as compared to the OA/ER algorithm. The disadvantage of both branch and bound and GBD is that these methods typically require the solution of many NLP subproblems during the search for the optimal solution (see Duran and Grossmann, 1986a). The OA/ER algorithm was developed to minimize the number of NLP subproblems to be solved since this is often the bottleneck in MINLP problems. Among other methods that are currently under development for solving MINLP problems is the feasibility technique by Mawengkang and Murtagh (1986). DICOPT As mentioned in the previous section, the OA/ER algorithm involves the alternating solution of NLP subproblems and MILP master problems. This algorithm has been implemented in the computer code DICOPT. MINOS (version 5.0 or 5.1) is used to solve the NLP subproblems and MPSX or ZOOM/XMP (Marsten, 1986) is used to solve the MILP master problems. A simplified flowchart of DICOPT is shown in Figure 1. First the algebraic formulation of the MINLP problem is supplied to GAMS which creates the formulation in a format recognizable to the NLP code (e.g. MINOS) and generates analytical gradients. The user must supply the initial values for the binary variables y and also a starting point for the continuous variables x which MINOS uses in solving the subproblem (NLP). The NLP solver is then called and returns a solution having one of three status conditions, optimal, infeasible, or unbounded* If the NLP subproblem is unbounded then the original MINLP problem is also unbounded and the OA/ER algorithm terminates. An infeasible status means that either the NLP subproblem is infeasible or that due to the starting point supplied to the solver, the problem appeared to be infeasible. Recall that the formulation (SNLP) can be used to avoid infeasibility due to inequality constraints and thus produce a point which satisfies the equality constraints while minimizing the violation of the inequalities. This is important when applying the equality-relaxation procedure because lagrange multipliers are required at a point which satisfies the nonlinear equalities in order to relax these equations to inequalities (see Kocis and Grossmann, 1987a). The third possibility is that the optimal (or at least locally optimal) solution to subproblem (NLP) was found. In this case, the upper bound on the solution to (MINLP) is updated as the minimum of the current upper bound and the objective function value from the NLP subproblem. Following the NLP subproblem a call is made to a program which formulates the MILP master problem. The master problem contains all linear constraints from (MINLP), an accumulation of linear approximations at previous iterations of the nonlinear constraints in (MINLP), and a set of integer cuts which eliminate from consideration the previously analyzed points y, k=l,2,..K. If the NLP subproblem solution is optimal or feasible with respect to the nonlinear equality constraints then first-order linearizations derived at x,y are included in the master problem. (Since gradient information is required by the NLP solver, this information is available for deriving the linearizations.) However, if for a given y the resulting NLP subproblem has no point x which satisfies the nonlinear equations, then the linear approximations are not added to the MILP formulation. The introduction of the integer cut will insure that the solution to the next master problem y + 1 is different from y. Having formulated the master problem, the next task in DICOPT is to call the appropriate MILP package. The MILP solver returns either the optimal integer solution (or an integer solution satisfying the specified tolerances) or an infeasible solution. If the status returned by the MILP solver is the optimal integer solution, the objective function value then provides a lower bound on the solution to (MINLP). If the master problem is infeasible or the predicted lower bound is greater than the current upper bound, then the OA/ER algorithm terminates and the optimal solution to (MINLP) has the objective function value of the current upper bound and the optimal solution point given by the corresponding NLP subproblem solution x* , y*. For the case when the lower bound from the master problem is less than the current upper bound, the binary variables y are temporarily fixed at the optimal value y from the master problem defining the next NLP subproblem. The optimal values of the continuous variables x from the master problem are used as the starting point for the next NLP subproblem meaning that the user is required to supply a starting point for only the first NLP subproblem. Since the linearizations are accumulated as iterations proceed, the MILP master problem provides an increasingly good approximation to the original MINLP. Hence, the quality of the starting points supplied through the master problem increases making the subsequent NLP subproblems easier to solve. Thus, the MILP master problem in the OA/ER algorithm not only predicts strong lower bounds, but it also provides excellent initial guesses for the NLP subproblems. The main components in DICOPT include GAMS, an NLP solver, an MILP solver, the code which formulates the master problem, and a control program which determines the flow between these four components. Efficiency and portability have been built into the computer implementation of DICOPT. Sparse matrix techniques as well as dynamic memory capabilities in GAMS are used to generate the linearizations in the MILP master problem. Currently, versions of DICOPT exist for IBM mainframes (CMS), VAX (VMS), and IBM-PC (DOS). Since MINOS and ZOOM/XMP are written in FORTRAN, they are available on these three systems. However, MPSX is available only for IBM mainframes. The code which formulates the master problem and the code which performs the necessary logic were written in FORTRAN. The command files, which are very small, are written in REXX for IBM/CMS, DCL for VAX/VMS, and BAT for IBM-PC/DOS. It should be noted that due to the robustness and efficiency of MPSX, the IBM version is the most suitable for large scale MINLP problems. To illustrate the use of the modelling language GAMS in DICOPT the algebraic formulation of a batch process design problem is given in Appendix C. Note the use of indexed sets which allows for a very compact problem formulation. In the solve statement given in the last line of the Appendix C listing, the command MIDNLP stands for the use of the OA/ER algorithm for solving the MINLP problem. Also, note that in DICOPT the user need not be concerned with the details of the MINLP algorithm (e.g. supplying the master problem formulation) as is the case of APROS developed by Paules and Floudas (1987). Here all the steps of the particular MINLP algorithm (OA/ER or GBD) must be supplied in the GAMS input file. COMPUTATIONAL RESULTS DICOPT, as implemented on the IBM/CMS system, has been tested on fifteen MINLP problems. A summary of problem size and computational results on an DBM-3090/600 (Cornell Theory Center) is given in Table L Problem REAC corresponds to the selection of a two reactor configuration, while problems PLAN and EX3 are small planning examples where discrete choices are made pertaining to candidate processes for producing a set of desired final products, CAP1 and CAP2 are examples of the pure-integer quadratic programming formulation of the capital budgeting problem (Kettani and Oral, 1987). Problem EX4 addresses the optimal positioning of a new product in a multiattribute space (Duran and Grossmann, 1986a). REL is a reliability optimization problem from Henley and Kumamoto (1985) and FLEX is an MINLP formulation of the flexibility problem of a heat exchanger network through the active set strategy of Grossmann and Floudas (1987). BATCH5 through BATCH12 correspond to optimal design problems for multiproduct batch plants (see Kocis and Grossmann, 1987b) and problem TFY involves the retrofit design of heat exchanger networks. Finally, the problem FLOW is an example of a structural flowsheet optimization problem (see Kocis and Grossmann, 1987b). The detailed formulations are given in Kocis (1988) and are available upon request from the authors. As shown in Table I, the size of the MINLP problems ranges from 2 to 60 0-1 variables, 0 to 410 continuous variables, and 1 to 421 constraints. Note that problem FLOW involves 140 nonlinear equations and problem EX4 25 nonlinear inequalities. The remaining problems exhibit fewer nonlinearities in the objective function or constraints. As can be seen, all the problems, except for the last one, required less than 30 seconds of CPU-time for solving both the NLP subproblems and MILP master problems. The time for the total overhead was always substantially less that the CPU-time for optimization. Also, the results in this table indicate that between 2 and 4 iterations are typically required. From the results in Table I, it is clear that the performance of the OA/ER in DICOPT algorithm is quite impressive considering both the size and complexity of the different test problems. Note that although the MILP master problems often require more CPU-time than the NLP subproblems (especially when few nonlinearities are present), they have the important role of not only reducing the number of iterations, but also to make the solution of subsequent NLP subproblems easier to solve by supplying good initial points. In fact, in all these problems, once the first NLP problem converged there was no difficulty in subsequent iterations. Finally, the results show that in problem FLOW, which involves complex NLP subproblems due to the large number of nonlinearities, the MILP master problems require much less time that the NLP subproblems. Another important point related to computational experience is the numerical solution of the MILP master problems. Since the constraints in the master problem are derived from linearizations, there is a tendency for the magnitude of the coefficients to cover a wide range of values making the numerical stability of the MILP solver important. Computational experience has indicated that MPSX outperforms ZOOM/XMP in solving the master problems both in efficiency and reliability. Results in Table II indicate that MPSX is more suitable for solving the MILP master problems as the size of the MINLP problem increases. In addition to the problems reported in Table I, the OA/ER algorithm has been applied to the design of gas pipelines (Duran and Grossmann, 1986b), retrofit design of batch processes (Vaselenak et al, 1987), and the structural flowsheet optimization with heat integration (Kocis and Grossmann, 1987a). In all these problems the OA/ER algorithm required between 2 and 5 iterations. Effective use of the OA/ER algorithm in chemical engineering applications has also been recently confirmed by other researchers. Floudas and Paules (1987) developed an MINLP formulation for the synthesis of heat integrated distillation sequences. With the OA/ER algorithm they were able to solve the resulting convex MINLP problem in 3 to 6 iterations. Harsh (1987) developed a mixed-integer programming approach to the retrofit design of existing plants. He successfully implemented the OA algorithm in the commercial process simulator FLOWTRAN and performed a retrofit design study of an ammonia plant in only 2 iterations. EXTENSIONS The OA/ER algorithm is guaranteed to find the global optimum of MINLP problems that involve convex objective function and inequalities, and quasiconvex relaxed equations. A difficulty, however, that may be encountered in MINLP problems is the presence of nonconvexities. In the OA/ER algorithm, nonconvexities can cause problems in two ways. First, since the procedure involves the solution of NLP subproblems, the presence of nonconvexities leads to the potential for local solutions in the subproblem (NLP). The global optimization of nonconvex NLP problems is an difficult task and rigorous algorithms currently do not exist The second difficulty that nonconvexities can lead to occurs in the MILP master problem of the OA/ER algorithm. The rigorous guarantee of the global optimum to (MINLP) is based on the master problem providing a valid lower bound on the optimal objective function value. However, when nonconvexities are present, the linearizations in the master problem may not necessarily provide valid outer approximations to the nonlinear feasible region. Therefore, the lower bounds predicted by the MILP master may not be valid lower bounds and the global optimum can be cut off (see Kocis and Grossmann, 1987b). In some special problems, nonconvexities can eliminated through convexifying transformations (e.g. log transformations as in Duran and Grossmann, 1986b). However, the identification of nonconvexities is often nontrivial for complex models and even if the nonconvexities are detected it is not always possible to convexity the problem. Thus, there is a need for a procedure which addresses the difficulty that nonconvexities cause in the MILP master problem. The important point to take note of is that convexity in the nonlinear functions is a sufficient condition to guarantee the global optimum. When these conditions are not satisfied the OA/ER algorithm may or may not obtain the global optimum solution. However, computational experience by the authors has shown that in a good number of cases the OA/ER algorithm will find the global optimum solution for nonconvex MINLP problems. This would then suggest that a suitable strategy to tackle these problems is to solve them in two phases as has been suggested recently by Kocis and Grossmann (1987b). In the first phase the OA/ER algorithm is applied in its original form but with special local and global tests for the identification of nonconvex functions that may cause the global solution to be missed. If nonconvexities are not detected the